R Notebook

data load

library(dplyr)

Attaching package: ‘dplyr’

The following objects are masked from ‘package:stats’:

    filter, lag

The following objects are masked from ‘package:base’:

    intersect, setdiff, setequal, union

library(data.table)

Registered S3 method overwritten by 'data.table':
  method           from
  print.data.table     
data.table 1.14.2 using 1 threads (see ?getDTthreads).  Latest news: r-datatable.com
**********
This installation of data.table has not detected OpenMP support. It should still work but in single-threaded mode.
This is a Mac. Please read https://mac.r-project.org/openmp/. Please engage with Apple and ask them for support. Check r-datatable.com for updates, and our Mac instructions here: https://github.com/Rdatatable/data.table/wiki/Installation. After several years of many reports of installation problems on Mac, it's time to gingerly point out that there have been no similar problems on Windows or Linux.
**********

Attaching package: ‘data.table’

The following objects are masked from ‘package:dplyr’:

    between, first, last

library(ggplot2)
library(ggmap)

Google's Terms of Service: https://cloud.google.com/maps-platform/terms/.
Please cite ggmap if you use it! See citation("ggmap") for details.

library(RColorBrewer)
library(leaflet)

Registered S3 method overwritten by 'htmlwidgets':
  method           from         
  print.htmlwidget tools:rstudio

library(scales)
library(plotly)

Attaching package: ‘plotly’

The following object is masked from ‘package:ggmap’:

    wind

The following object is masked from ‘package:ggplot2’:

    last_plot

The following object is masked from ‘package:stats’:

    filter

The following object is masked from ‘package:graphics’:

    layout

library(rlist)
library(broom)

powerplants=read.csv('global_power_plant_database.csv')

powerplants[is.na(powerplants)] = 0
powerplantsCountry = powerplants %>% filter(country_long == "United States of America")

Data prep

powerplantsCountryFuel = powerplantsCountry %>% group_by(primary_fuel) 

powerplantsCountryFuelYear = powerplantsCountryFuel %>% summarise(energy2013 = 0.5* (sum(generation_gwh_2013) + sum(estimated_generation_gwh_2013)),
          energy2014 = 0.5* (sum(generation_gwh_2014) + sum(estimated_generation_gwh_2014)),
          energy2015 = 0.5* (sum(generation_gwh_2015) + sum(estimated_generation_gwh_2015)),
          energy2016 = 0.5* (sum(generation_gwh_2016) + sum(estimated_generation_gwh_2016)),
          energy2017 = 0.5* (sum(generation_gwh_2017) + sum(estimated_generation_gwh_2017)),
          )

powerplantsCountryFuelYearDF = as.data.frame(powerplantsCountryFuelYear)
powerplantsCountryFuelYearLabeled = powerplantsCountryFuelYearDF[,-1]
row.names(powerplantsCountryFuelYearLabeled) = powerplantsCountryFuelYearDF[,1]
transposePowerplantsCountryFuelYear = transpose(powerplantsCountryFuelYearLabeled)
colnames(transposePowerplantsCountryFuelYear) = powerplantsCountryFuelYearDF[,1]
transposePowerplantsCountryFuelYear$Year = c(2013,2014,2015,2016,2017)


totalEnergy = powerplantsCountryFuelYear %>% mutate(total = energy2013+ energy2014 +energy2015 + energy2016 + energy2017)

map

powertype = "Wind"
leaflet(powerplants %>% filter(primary_fuel == powertype)) %>% addTiles() %>% addMarkers(~longitude,~latitude,popup= ~primary_fuel,label= ~name)

Total energy bar chart

totalEnergyChart=ggplot(totalEnergy,aes(x=reorder(primary_fuel,total),y=total,fill=primary_fuel)) + geom_bar(stat="identity") + coord_flip() + scale_y_continuous(trans= 'log10', labels = trans_format('log10',math_format(10^.x))) + theme(axis.text.x = element_text(angle=90)) + ylab("Energy Accumulated MWH") + xlab("Primary Fuel")
totalEnergyChart

energy produced

data_long = melt(transposePowerplantsCountryFuelYear, id="Year")

Warning in melt(transposePowerplantsCountryFuelYear, id = "Year") :
  The melt generic in data.table has been passed a data.frame and will attempt to redirect to the relevant reshape2 method; please note that reshape2 is deprecated, and this redirection is now deprecated as well. To continue using melt methods from reshape2 while both libraries are attached, e.g. melt.list, you can prepend the namespace like reshape2::melt(transposePowerplantsCountryFuelYear). In the next version, this warning will become an error.

colnames(data_long) = c("Year","primary_fuel","energy")

lineplot = ggplot(data_long, aes(x=Year,y=energy, color= primary_fuel)) + ylab("Energy Produced MWH") + xlab("Year") + geom_line() + scale_y_continuous(trans= 'log10', labels = trans_format('log10',math_format(10^.x))) +ggtitle("Energy Produced over the Years per Type")

ggplotly(lineplot)

Warning in is.na(ticktext) :
  is.na() applied to non-(list or vector) of type 'expression'

number of plants by type

numberOfPlants = powerplantsCountryFuel %>% summarize(number = n())
ggplot(numberOfPlants,aes(x=reorder(primary_fuel,number), y=number,fill=primary_fuel)) + geom_bar(stat="identity") + theme(axis.text.x = element_text(angle=90)) + coord_flip() + scale_y_continuous(trans= 'log10', labels = trans_format('log10',math_format(10^.x))) + xlab("Primary Fuel") + ylab("Number of Plants") +ggtitle("Number of Plants per Type")

powerplantsCountryFuelAndNumber = merge(totalEnergy,numberOfPlants,"primary_fuel","primary_fuel")
powerplantsCountryFuelAndNumber = powerplantsCountryFuelAndNumber %>% mutate(EnergyPerPlant = total/number)

totalEnergyPerPlantChart=ggplot(powerplantsCountryFuelAndNumber,aes(x=reorder(primary_fuel,EnergyPerPlant),y=EnergyPerPlant,fill=primary_fuel)) + geom_bar(stat="identity") + coord_flip() + scale_y_continuous(trans= 'log10', labels = trans_format('log10',math_format(10^.x))) + theme(axis.text.x = element_text(angle=90)) + ylab("Average Energy Accumulated MWH per Plant") + xlab("Primary Fuel")
totalEnergyPerPlantChart

powerplantsCountryFuelTotal = powerplantsCountryFuel %>%
  mutate(totalEnergy = 0.5 * (
    generation_gwh_2013 + estimated_generation_gwh_2013 +
    generation_gwh_2014 + estimated_generation_gwh_2014 +
    generation_gwh_2015 + estimated_generation_gwh_2015 +
    generation_gwh_2016 + estimated_generation_gwh_2016 +
    generation_gwh_2017 + estimated_generation_gwh_2017
    ))
totalEnergyBoxPlot = ggplot(powerplantsCountryFuelTotal,aes(x=reorder(primary_fuel,totalEnergy),y=totalEnergy,fill=primary_fuel)) + geom_boxplot(show.legend = FALSE) + coord_flip() + scale_y_continuous(trans= 'log10', labels = trans_format('log10',math_format(10^.x))) + theme(axis.text.x = element_text(angle=90)) + ylab("Average Energy Accumulated MWH per Plant") + xlab("Primary Fuel")
totalEnergyBoxPlot

Warning in self$trans$transform(x) : NaNs produced
Warning: Transformation introduced infinite values in continuous y-axis
Warning: Removed 117 rows containing non-finite values (stat_boxplot).

windPowerplants = powerplantsCountryFuelTotal %>% filter(primary_fuel == "Wind")
oilPowerplants = powerplantsCountryFuelTotal %>% filter(primary_fuel == "Oil")

# Compute the analysis of variance
res.aov <- aov(totalEnergy ~ primary_fuel, data = powerplantsCountryFuelTotal)
# Summary of the analysis
summary(res.aov)

               Df    Sum Sq   Mean Sq F value Pr(>F)    
primary_fuel   13 1.253e+11 9.641e+09     762 <2e-16 ***
Residuals    9819 1.242e+11 1.265e+07                   
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

thsd = TukeyHSD(res.aov, "primary_fuel", ordered = TRUE)
tidy(thsd) %>% filter(contrast=="Wind-Oil")